Matter-wave solitons supported by dissipation 
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We show how novel types of long-lived self-localized matter waves can be constructed with Bose- 
Einstein condensates. The ingredients leading to such structures are a spatial phase generating a 
flux of atoms towards the condensate center and the dissipative mechanism provided by the inelastic 
three-body collisions in atomic Bose-Einstein condensates. The outcome is an striking example of 
nonlinear structure supported by dissipation. 
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The realization of Bose-Einstein condensation with ul- 
tracold atomic gases has opened the door for many spec- 
tacular realizations of matter waves. In particular, the 
elastic two-body collisions between condensed atoms pro- 
vide an effective nonlinear interaction within the atomic 
cloud. These nonlinear interactions have been used to 
obtain different types of self-localized matter waves Q 
denoted as matter- wave solitons. 

Most proposed and all experimentally realized solitons 
in ultracold atomic gases lead to self-localization along 
one spatial dimension while in the others they must be 
externally confined by magnetic or optical means. Truly 
multidimensional solitons cannot be supported by the 
nonlinear interactions arising from elastic two-body col- 
lisions and nonlinearity is not able to balance robustly 
dispersion these scenarios H Q . 

In Bose-Einstein condensates (BECs) with Feschbach 
resonance management it has been predicted that cer- 
tain two-dimensional solitons may exist 0, H Q but the 
idea does not extend trivially to three spatial dimensions 
without the adition of extra (external) confining poten- 
tials. 

Another physical effect coming from the interaction 
between atoms in the condensate is dissipation. Three- 
body collisions usually lead expulsion of atoms from the 
condensate leading to an effective dissipation which is 
thought to be responsible for its finite lifetime. 

Dissipation usually acts against self-localization since 
it tends to take the system closer to the linear situation 
where no stable nonlinear structures exist. However, in 
some contexts dissipation has been shown to play an sta- 
bilizing role. In Ref. pj the addition of a phenomenolog- 
ical Landau damping to the Feschbach-resonance man- 
aged model has been shown to enhance stability. Another 
example is nonspreading (linear) wave packets with ex- 
ternal imaginary potential Q. 

In this letter we present novel nonlinear structures in 
Bose-Einstein condensates self-trapped by the effect of 
nonlinear dissipative terms. The physical idea behind 
this paper is that imprinting an appropriate spatially de- 
pendent phase on a BEC leads to a flux of particles from 
its periphery to the center which compensates the parti- 



cles lost by three-body inelastic collisions. The outcome 
is a long-lived nonlinear structure supported by dissipa- 
tion. 

We will work in the mean field approximation in which 
a BEC is modelled by the Gross-Pitaevskii equation 

ot 2m m 
where V accounts for any external trapping potential, 
a s is the s— wave scattering length for elastic two-body 
collisions and T3 = — HK^/ 12 , K3 being the threshold- 
less three-body recombination rate 9] . Eq. (QJ has been 
used as a model to study collapse in BECs with attrac- 
tive interactions El El El El 111 El El However, the 
nonlinear structures to be studied in this paper are not re- 
lated with collapse phenomena, a fact that we will stress 
by working out examples with a s = or even a s > 0. 
Models similar to Eq. Q arise in the propagation of op- 
tical beams in Kerr media with multiphoton absorption 
processes [HtHI- 

First, we change to the dimensionless variables r = 
x/clq, t = t/T (z^o = 1/T) and i/j = \/^o^ where = 
ft/a^ra, and ao is a characteristic size of the BEC so that 
Eq. Q becomes (from now on we set V = 0) 
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with g = 4ira 3 /ao and 7 = —K3/ {12vq<1q). 

To obtain stationary self-trapped solutions of Eq. (0 
we first write it in the modulus-phase representation, i. 
e. ip(r,T) = A(r, t) exp[i<p(r, r)], with A(r,r) > 
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Stationary solutions of Eq. @ satisfy d T A 2 = and 
from Eq. (JHh) we get <p(r, t) = (fi(r) — St, with 5 > 0. 
Rescaling the spatial variables with p = \f2br we find 
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In this paper we will focus on solutions with radially sym- 
metric amplitude A(p) = R(p). In quasi- two dimensional 
situations we will study solutions with phase given by 
4>{p,t) = *&(p) +m0, being = atan(p2/pi) i- e - includ- 
ing a possible vorticity. In three spatial dimensions we 
will restrict our attention to spherically symmetric sta- 
tionary solutions 4>{p,t) = Then, Eq. (|4a[) can be 
integrated by using the divergence theorem and we get 
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R 2 p d ~ 1 S 

Eq. © allows us to obtain the asymptotic behavior 
of the phase, &'(p) ~ j^R 4 (p)p J when p <C 1, and 
~ Q/ [/ _1 ^ 2 (P)] 5 f o r P > 1, where we assume 
Q = j Jq 00 rd ~ X R Q dr to be a finite quantity. The ampli- 
tude equation (|4b| becomes 
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R - tR 3 + ii = 0, 



(6) 



plus i?(oo) = and R(0) = R (m = 0) or R'{p) 
(m 7^ 0). In the later case, it is easy to obtain R(p) ~ 
i?'(0)p |m| for small p. For large p we get 



R" + d=l& + qR~0, 



(7) 



which gives R(p) ~ asymptotically. 

To obtain the following order in the approximation of 
R(p) for large p we define = R{p)p^ d ~ 1 ^ 2 , which sat- 
isfies the Ermakov-Pinney equation h" + h = Q 2 /h 3 
whose solution can be found in explicit form 



R(p) 



^ cos 2 p + (Ci cos p + Ci sin p) 2 



p (d-l)/2 



p » 1. (8) 



When d = 2 this result is similar to that of Ref. ^| 
with the different definition of Q. Eq. JHJ) guarantees the 
finiteness of Q which was assumed previously. 

We have solved numerically Eqs. (|la|) and © to find 
specific stationary configurations. Some results for g = 
are shown in Fig. for m = [Figs. GJa) and [2] and 
m = 1 [Fig. IHb)] as a function of the shooting param- 
eter a = 5/7i? 4 (0). Solutions are found above a cer- 
tain threshold value a > a t h = 0.54, 2.93, 0.36 corre- 
sponding to the 2D solution without and with vorticity 
(m = 1), and 3D solution, respectively. When the 
effect of the nonlinearity is included (see Fig. OJ) we ob- 
tain a compression of the solutions for g < [Fig. Efa)] 
and an expansion leading to slower decay for g > [Fig. 
Etb)]. The threshold a t h is lowered when g < Oand in- 
creased when g > there being a maximum positive 
value g* — 5/R 2 (0) above which stationary solutions do 
not exist. It is remarkable that stationary solutions sup- 
ported by dissipation can exist even when the scattering 
length is positive. In three spatial dimensions we find 
a similar behavior with a faster decay of the amplitude 
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FIG. 1: Radial profile of stationary solutions of Eq. @ with 
g = and (a) a = 8/[<yR 4 (0)] = 0.55, 0.65, 0.8, oo (from 
higher to lower ones), (b) Stationary radial profile of vortex 
solutions with m = ±1 for 5/j[R'(0)] 4 = 2.9, 3.2, 4.1, oo. 
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FIG. 2: Phase and its derivative for solutions of Eq. @ with 
g = and (a-b) a = 5/ [jR 4 (0)] = 0.65 (c-d) a = 5. 
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FIG. 3: [Color online] Radial profile of solutions of Eq. @ 
with (a) S/[jR 4 (0)] = 0.8 and gR 2 (0)/6 = -1.3 (dotted 
line), (dashed line), and 0.35 (solid line), (b) 5/ [yR 4 (0)] = 
5, gR 2 (0)/S = -2 (dotted line), (dashed line), and 0.85 
(solid line). 




FIG. 4: [Color online] Radial profile of solutions of Eq. @ 
for d = 3 with (a) 5/ [jR 4 (0)] = 0.36, 0.4, 0.5, oo and g = 
(b) 5/ [^R 4 (0)] = 0.4 and gR 2 (0)/S = -2 (dotted line), 
(dashed line), and 0.1 (solid line). 
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FIG. 5: [Color online] Time evolution of a stationary solution 
of Eq. Q for 7 = -2 x 10" 10 (from K 3 = 2 x 10" 26 cm 6 /s 
for 7 Li), g = 0, S = 0.008, R = 77. The solution is solved 
in a square region of 1.2 mm size and the stationary solution 
is set to zero outside a disk of radius 450 /xm which leaves 
about No = 5.7 x 10 7 atoms inside, (a-c) Pseudocolor plot 
showing \ip(x,y,t)\ for different times, (d-e) Evolution of the 
(d) amplitude of the wavefunction and (e) number of particles 
in the condensate. 



(Fig. [3J. In all cases the oscillatory behavior of the so- 
lutions and their decay are well described by Eq. (JSJ 

How can stationary solutions exist in a system with dis- 
sipation? From the previous asymptotic analysis R(p) ~ 
l/p( d - 1 )/ 2 j thus we can estimate the number of particles 
as N ~ J °° r d ~ 1 R 2 dr = 00. Moreover 



= 2 7 f |Vfdr < 00. 
at J Rd 



(9) 



Intuitively speaking, it is the fact that the solution has 
an infinite number of particles together with a finite par- 
ticle loss given by the right hand side of Eq. (0 which 
allows for the existence of stationary solutions even in 
the presence of three-body recombination losses. 

However, realistic solutions must have a finite number 
of particles. A simple way to construct realistic quasi- 
stationary solutions is to cut the stationary solutions at 
an specific point far enough from the origin. In Fig. Owe 
show simulations of the evolution of one of those modified 
stationary solutions. The central soliton survives for very 
long times with quasi-stationary amplitude as shown in 
Fig. EJd). This is a consequence of the phase structure 
of stationary solutions (see Fig. EJ), which leads to a par- 
ticle flux towards the condensate center and thus induce 
a refilling mechanism of the central soliton which is man- 
ifest in Fig. (a-c). The rings surrounding the central 




FIG. 6: [Color online] Evolution of a stationary solution of 
Eq. (0) computed with 5 = 10~ 2 , 7 = -2 x 10~ 10 . Shown are 
the evolutions of the amplitude (a,c) and number of particles 
(b,d) in the condensate for (a-b) g = —5.7 x 10 -6 , Ro = 132, 
on a square of side 800 /xra, with the solution cut at 300 jLtm, 
leading to iV = 7.6 x 10 7 atoms and (c-d) g = 7.4 x 10~ 7 , 
Ro =82, on a square physical domain of side 1.2 mm with 
the profile set to zero outside a disk of radius 450 /xm which 
implies No = 1.2 X 10 8 atoms. 



peak play the role of a reservoir of atoms and constantly 
feed the central peak and dissapear progressively as time 
proceeds. From the practical point of view and since the 
amplitude of the rings is small what one would observe 
is a very long lived soliton lasting for times of the order 
of the condensate lifetime. 

The number of atoms used in our simulations are about 
an order of magnitude above those present in realistic 
experiments, however what it is really important is the 
shooting parameter 5/ [7i? 4 (0)j. By choosing smaller 
values for R(0) one may get smaller particle numbers 
but at the price of decreasing S thus changing the solu- 
tion and widening the spatial scales. These are options 
which must be taken on the basis of particular experi- 
mental scenarios. 

We have checked that the same mechanism works for 
vortex type solutions as shown in Fig. [7| 

The scenario described above is also valid for nonzero 
(positive and negative) values of g as shown in Fig. [SJ 
This means that the phenomenon described here is not 
a result of elastic two-body collisions. Moreover it seems 
that the phenomenon results is a very long-lived cen- 
tral soliton even in condensates with positive scattering 
length, which is a very surprising result. 

We have not found any signs of the azimuthal instabil- 
ities predicted in Ref. 18] for systems with higher order 
dissipation. This fact requires further investigacion and 
could be either because of a suppression of this instability 
in our situation or just a consequence of a finite lifetime 
of our structure, smaller than the typical time for the 
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FIG. 7: [Color online] (a-c) Pseudocolor plots showing the 
evolution of \ip(x, y, t)\ 2 for different time values of a vortex 
type,ra — 1, solution computed with g — 0, 5 = 10 -2 , Rmax — 
57, on a square of side 900 /im, and setting the profile to zero 
outside a disk of radius 350 /im. The total number of atoms 
is No = 7.6 x 10 7 . 

instability to set in. 

Experimental generation. - To generate these structures 
in real experiments it would be very important to obtain 
appropriate initial data with amplitude and, more im- 
portantly, phase close to those corresponding to station- 
ary solutions. In quasi- 2D condensates phase imprinting 
methods [2(1 0, can be used. In fact, imprinting 

radially symmetric phases is simpler to do using absorp- 
tion plates than originally proposed for vortices . As 
to the amplitude profile it could be achieved by using a 
Bessel beam instead of a gaussian beam for trapping the 
condensate before releasing it. 

In fully three-dimensional Bose-Einstein condensates 
phase imprinting methods would be more difficult to ap- 
ply. However, three dimensional systems with attractive 
scattering length and three-body recombination in the 
collapsing regime develop the phenomenon known as su- 
perstrong collapse 0, 0> EH which has been ignored in 
the physical literature studying collapse in Bose-Einstein 
condensates. This means that a collapsing condensate 
spontaneously develops a structure of shells in both the 
amplitude and phase similar to those shown in Fig. E] 
constantly feeding atoms to the region of higher density 
of atoms. In fact, those structures are present in numer- 
ical simulations of collapse in Bose-Einstein condensates 
[HI EM E3 • The phenomenon of spontaneous phase self- 
modulation in superstrong collapse provides a way for 
generating the matter- wave solitons presented in this pa- 
per. If during the initial stage of a blow-up event the 
scattering length is changed to a subcritical value the vi- 
olent compression and high losses associated to collapse 
would be suppressed but the superimposed phase struc- 
ture could evolve to one of the attractors of the system, 
i.e. the stable stationary states. This procedure could be 
a way to generate the matter wave soliton supported by 
dissipation proposed in this paper. 

This work has been supported by grants BFM2003- 
02832 (Ministerio de Educacion y Ciencia, MEC) and 
PAI-05-001 (Consejeria de Educacion y Ciencia de la 
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